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Abstract 

We investigate the quantum behavior of a mesoscopic two-boson sys- 
tem produced by number-squeezing ultracold gases of alkali metal atoms. 
The quantum Poincare maps of the wavefunctions are affected by chaos 
in those regions of the phase space where the classical dynamics produces 
features that are comparable to h. We also investigate the possibility for 
quantum control in the dynamics of excitations in these systems. Con- 
trolled excitations are mediated by pulsed signals that cause Stimulated 
Raman Adiabatic passage (STIRAP) from the ground state to a state of 
higher energy. The dynamics of this transition is affected by chaos caused 
by the pulses in certain regions of the phase space. A transition to chaos 
can thus provide a method of controlling STIRAP. 



1 Introduction 

There have been significant advancements in techniques for cooling and trapping 
ultracold atom gases in recent years, facilitating the experimental realization of 
Bose-Einstein condensation in dilute alkali gases (specifically 87 Kb and 23 Na) 
in 1995 pQ [5] [5] [I]- Since then, numerous studies of these condensates have 
been accomplished. In addition, experiments have been conducted that have 
obtained boson systems in number squeezed states from ultracold gases of alkali 
atoms in optical traps [S]. Thus, it is now possible to create a mesoscopic two- 
boson system of Sodium or Rubidium. Individual alkali atoms from a BEC 
reservoir can be subjected to quantum tweezers (Gaussian lasers exploiting the 
Landau-Zener tunnelling between the reservoir levels and the levels in the laser 
beam) |6 . Other recent methods include number-squeezing the BEC itself by 
"culling" atoms from a trapped condensate down to a sub-poissonian regime, 
making the number uncertainty small enough to be ignored [5]. 

A two-boson system can be subjected to a micrometer-scale double well po- 
tential using coherent laser beams. The methods that are applicable in the 



required length scales are overlapping cross-sections of Gaussian lasers [7] . Nu- 
merous methods can be used to control the extent of the overlap, such as an 
Acousto-Optical Modulator (AOM) vibrating with two sound waves [3]. Each 
sound wave causes an impinging laser beam with a Gaussian cross section to 
Bragg-diffract at an angle. The separation between the two diffracted beams can 
be adjusted with relative ease, and consequently so can the length scale of the 
double well that is generated by focussing the two diffracting beams into parallel 
beams (the double well lies in the direction lateral to the propagation) [S] . Other 
applicable techniques for double well generation are small volume optical traps 
formed by multiple lasers propagating as Gaussian sheets [5], and blue-detuned 
far-off-resonant laser light to add a potential hill in the middle of a cigar-shaped 
magnetic trap [ID] . 

A two-boson system, in an optical trap, can be subjected to stimulated 
Raman scattering. As we shall show, coherent population transfer from the 
ground state into one of the excited states can be achieved using radiation 
pulses. Time-modulated (i.e pulsed) electric fields or laser radiation pulses, 
applied sequentially in a counter intuitive manner, can be used to facilitate 
this process. If the time scale of the pulse modulation is sufficiently large, the 
Raman process is adiabatic (called Stimulated Raman Adiabatic Passage or 
STIRAP) [TT] [12]. 

The analysis of traditional STIRAP systems involves relatively weak radia- 
tion pulses and the rotating wave approximation on the three involved levels [12] . 
More recently, multilevel transitions have been performed and Floquet theory 
has been used to tell us how the population evolves in time [13] [H] [15] [16] , 

We investigate the appearence of avoided crossings contributed by resonat- 
ing levels other than the ones that the radiation pulses connect. These avoided 
crossings in the Floquet eigenphases appear due to level repulsion caused by a 
loss of symmetry/degeneracy (actual crossings) [17], and affect the statistical 
properties of the spectra, bringing them close to that predicted by random ma- 
trix theory. These are connected with the dynamics of the underlying classical 
system, which undergo a transition from KAM tori to chaos in this region of the 
paramter stace |17j . Thus, this work also demonstrates the quantum effects of 
chaos, induced by the radiation, on multilevel transitions in a 2-boson system. 

In the following sections, we describe the behavior of two bosons in a one- 
dimensional double-well potential which we shall model as detailed in Section 
2. In Section 3 we will describe the classical dynamics, where the pseudopoten- 
tial interaction in one-dimension is approximated by a Gaussian potential. In 
Section 4, we will discuss the quantum eigenstates of this system and compare 
the quantum phase space of the eigenstates with the classical phase space. In 
Section 5, we will proceed to describe the dynamics of excitations of this system 
which will be driven by two sinusoidal electric field pulses applied in sequence. 
The frequencies of the fields are chosen to connect specific undriven energy 
levels of the system. Section 5.1 poses the Schrodinger equation for the quan- 
tum dynamics for this system. Section 5.2 introduces Floquet theory and the 
numerical methods implemented to evaluate the Floquet matrix. Section 6 dis- 
cusses a specific set of parameters for the general STIRAP dynamics discussed 
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in Section 5. In subsections 6.1-6.3, we look at the eigenvalue spectra of the 
Floquct matrices for two amplitudes and use them to explain the actual dynam- 
ics of the system obtained from ab initio numerical solutions of Schrodinger's 
equation with the same parameters. We repeat these numerical methods for 
a different set of parameters in Section 7, where the wells are slightly deeper 
and an additional resonance exists in the eigenvalue spectrum which affect the 
coherent excitations. For large pulse amplitudes, the dynamics of the STIRAP 
process is also noticeably affected by the presence of chaos in the stronger pulse 
amplitudes in both cases. Concluding remarks are made on Section 8. 



2 The Basic Model 

Our system consists of two alkali metal bosons confined to a double-well optical 
potential. The two bosons can be obtained from a cold-atom system confined in 
a Magneto-Optical Trap (MOT) by number squeezing, followed by laser culling 
in an optical dipole trap (ODT). If the optical laser is far-detuned from the 
internal atomic resonances, they can be treated as point particles. The effective 
interaction between the bosons, in three dimensions, is obtained in the long 
wavelength approximation to be 

3d/ \ 47r?i 2 a s . 

U dCl Xi-X2) = <5xi-x 2 , (1) 

m 

where Ti is Planck's constant, a s is the s-wave scattering length and x^ = 
(xi, yi, Zi) is the displacement of the ith particle [TS] [IS]. Therefore, the energy 
eigenvalues of the two particle system are given by the Schrodinger equation 



* Sj (x 1 ,x 2 ) 
+ [Vt(xi) + V T (x 2 )] 9 Bj (xi,x 2 ) 

m 

where m is the mass of the particles, V e (xi) is the optical double- well potential 
and depends only on a single coordinate, and Vt(x^) is the potential profile of 
the MOT. ^>e (xi,x 2 ) = (xi, x 2 |-Ej) is a symmetrized energy eigenstate of the 
interacting two-particle system with energy Ej . 

The system can be confined in two spacial (radial) directions so that the 
essential dynamics occurs in the x - direction. This can be achieved by using 
anisotropic magnetic traps with high aspect ratio [20], where the radial trap 
frequencies are considerably larger than the axial frequency |21) . In that case, 
the radial part of the wavefunction is is not affected by the interaction and the 
energy eigenstate can be decomposed into an axial part (along the x-axis) and 
a radial part such that 

* jEj .(xi,x 2 ) = tp J (x 1 ,x 2 )(j)±{xi± 7 X2±), (3) 
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where Xi^ = (yi,Zi) and 4 > x( x i±^ x 2x) denotes the noninteracting ground state 
of the system in the {y, z} plane. If we assume the trap to be a radial harmonic 
oscillator of stiffness ui s then, multiplying Eq. @ by (xu_ , X2j_) and integrat- 
ing in all radial coordinates reduces the Hamiltonian to a form which describes 
motion along the axial direction, 



2 2 

|i- + V e ( Xl ) + ^ + V e (x 2 ) + u( Xl 
2m Am 



X 2 ) 



t/)j(xi,X2) = Ejipj(xi,x 2 ), (4) 



where 



u(x\ — x 2 ) = Aa s uJshS(xi 



X 2 ) 



(5) 



and the energy Ej contains contributions from the radial ground state. Thus, 
the system is effectively one dimensional. 

We will consider the case of two identical bosons confined to a quartic double 
well potential. The total Hamiltonian for the system is 



Pi 



Pi 



H = -P- + -P- + v(-2ax\ + bx\) + v(-2ax\ + bxi) + 4,a s uj s h8{x 1 - x 2 ), (6) 
2m 2m 

where pi is the momentum of the ith particle (i— 1,2), Xi is the position of the 
ith particle along the x-axis, v determines the depth of the double well potential, 
and a and b determine its spatial extent. 

It is useful to write the Hamiltonian in terms of dimensionless parameters 
(p-, x^H'). We introduce a unit of length, L u — 50nm, which (as we will see) 
is appropriate to the systems we consider here. Then let x^ = jf^, a = L~ 2 , 



T ~ 4 n' 

7 Pi 



E u — 



and T„ 



2mLl 



- and t' = if where P u = j r , 



2mL'i 

parameters, we obtain 



E u ' ^ u — E. 

If we now drop the primes on the dimensionless 



H = pf +vl + V (-2xl + x\) + V (-2x 2 2 + x 4 2 ) + U 5(xi - x 2 ) 



(7) 



The unit of length, L u = bOnm, is two orders of magnitude below typical length 
scales in double wells created with Gaussian lasers [7] , as well as Gaussian sheet 
lasers, or lasers diffracted off of acousto-optical modulators vibrating with two 
sound waves [5]. Nonetheless, this length scale should be attainable in these 
setups. For instances, well dimensions in a double well generated using an AOM 
are proportional to the difference between the acoustic frequencies (typically 10 
MHz [9]), producing a barrier length in the micrometer range. It should be 
relatively straightforward to reduce that frequency difference by two orders of 
magnitude and adjusted it to produce a double well of a barrier length of lOOnm, 
or an L u of 50nm, as is required here. The lifetime of the magnetic traps in 
these experiments is approximately 3 — 4 s. If we use 87 Rb alkali atoms as our 
bosons then the value of T u is about 7/is, making typical trap lifetimes (about 
5s) translate to 800,000 units of T u . The characteristic energy scale E u is very 
small, and corresponds to photons with frequency of about 2AKHz. Therefore, 
monochromatic coherent radiation at around this frequency is necessary for 
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STIRAP excitations, and can be obtained using noble gas masers [8]. Fig. [T] 
shows a plot of the quartic double well V(x) = Vo(— 2x 2 + x A ) for well depth 
V a = 4.91345043. 

3 Classical Dynamics of the Undriven System 

In order to study the classical dynamics of the system governed by the Hamilto- 
nian in Eqn. j7]) the contact potential U(x\ — xq) = Uq5(x\ — X2) can be replaced 
by a Gaussian shaped potential of suitably chosen width a c such that 

U G {x x -x 2 ) = -^=e (8) 

(T C V Z7T 

We have noticed no discernable difference in the quantum case between using 
the Gaussian for the interaction and using the delta function provided that the 
width of the Gaussian a c is sufficiently small. Too small a width generates 
unresolvable errors in the numerical solution to the classical dynamics and so 
an optimum width was chosen at a c — 0.005. 

In the absence of a time dependent external field, the system is conservative, 
has two degrees of freedom, and energy Eq is constant. Therefore, the system 
is confined to a three dimensional surface in a four dimensional phase space. 
Fig. [2] shows Poincare surfaces of section of the classical phase space of this 
system for Vb = 4.91345043 and Uq = —1.0 (attractive interaction) and Eg 
fixed. The momentum pi and displacement x\ of particle-1 are plotted each 
time the trajectory of particle-2 crosses the point x-z = 1 with p2 > such that 
Eq is fixed and therefore they show the behavior of momentum and position of 
particle-1 for that energy. Energy conservation places bounds on the values of 
the coordinates p\,p2,Xi, xi so the trajectories are confined to a finite region of 
the phase space. 

3.1 Relevant Energies of the Quantum System 

Figs. [5]a through Oc are surfaces of section for three different energies E\ = 
—3.7195, E2 = —2.66655 and E4 — 2.5986, respectively. The energies were 
chosen to match quantum energy levels that will be connected by STIRAP 
pulses. The numerical integration was done by the 4 th order implicit Rungc 
Kutta (Prince Dormand) method [22] using the appropriate subroutine in the 
GNU Scientific Library 23J. The integration was done non-adaptively, with 
a fixed temporal stepsize of 10 -3 . In Fig. [2]a (E = —3.7195), we can see 
several regions of interest. There are three kinds of dynamics at this energy 
vis-a-vis the relative energies of each particle. In the case that the particles 
are in separate wells, trajectories exist where they don't see each other and are 
therefore the same as that of a single particle in a double well. They are seen 
as large KAM (Kolmogorov-Arnold-Moser) tori around x\ = —1 in Fig. HJa.l. 
These trajectories lie between two chaotic regions. Both chaotic regions are 
produced by the case when the energy of one particle is set to a positive value 
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(thus taking it above the wells), and the energy of the other particle is decreased 
so that they both add up to E (the particles being kept sufficiently far apart 
at t = so as to make the interaction negligible at that time). The resultant 
dynamics is chaotic due to the interaction experienced by the particles when 
they approach each other during motion. The cases when the particle being 
strobed has high momentum cause the chaos at the separatrix coupling both 
wells. The island seen immediately around x\ = —1 in Fig. [2]a.l is the result 
of a bifurcation that occurs at lower energies, and will be discussed in the next 
section. 

The case when both particles are in the same well are seen as the highly 
elongated tori around x\ = 1 in Fig. [2]a.2 due to the interaction between the 
particles. In Fig. [2jb.l, we note the disappearance of the bifurcated island 
immediately around x\ = —1, the result of a bifurcation in reverse (as energy is 
increased). The left-right asymmetry observed at lower energies is reduced at 
the energy increases and (apart from the elongated tori), disappears in Figs.Oc. 
The region of chaos seen deep inside the potential wells of Fig.[5]a.l merges with 
the region of chaos at the separatrix for higher the value of energy in Fig. [51c. 1. 

3.2 Lower Energies 

Figs. [3la through Of are surfaces of section for six energies, all less than the 
quantum ground state energy, and shown in decreasing order in energy. Only 
the x\ < half of the phase space is shown. For sufficiently low energies, the 
two particles either oscillate in two wells independently or together in the same 
well. The former case is seen in Fig.[3]f where the periodic motion of one particle 
about xi = —1 is visible at energy E = —5. As we increase the total energy, it 
becomes possible for one particle to break it's confinement for certain values of 
it's initial coordinates and influence the dynamics of the other particle through 
the interaction. This produces prominent chaotic behavior as seen in the figures. 
Increasing the total energy further from E = —4.5 to E — —4.0 (Figs.[3]d to[3]c) 
produces a bifurcation as the total energy is increased further (see Figs. [3jc to 
la). 

4 Quantum Mechanics of the Interacting Sys- 
tem 

In this section, we discuss the quantum mechanics of two interacting bosons 
whose Hamiltonian is given by Eq. (|7|). We first discuss the basis used to 
diagonalize the Hamiltonian. We then show configuration space and phase space 
plots of the key eigenstates of the system. 

4.1 Diagonalization of the Hamiltonian 

In order to diagonalize the Hamiltonian in Eq. ([7]) , we use the eigenstates of two 
noninteracting bosons in a hard-wall box as a representation to formulate the 
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matrix elements. We choose a box of width L = 3.5 (in dimensionless units) so 
that an adequate balance is achieved between truncability and accuracy. The 



Hamiltonian for a single particle in a box is h 
eigenstates are given by 



p y x < \L\ and the energy 



$n(x) 



1 



in, X) 



sm 



mr x ^ 
~2 [ L ~ 1 ' 



(9) 



where n = 1, 2, ...oo. 

The 2-particle boson states of the box system are obtained by symmetrizing 
the 2-particle states to obtain a complete orthonormal basis of symmetrized 
2-boson states: 



(xi,x 2 \ni,n 2 ) (s) 



1 



(a;i|m)(a;2|n2) + (xi\n 2 }(x 2 \ni)}. (10) 



These states are then used to create a Hamiltonian matrix from Eq. ([7|. The 
eigenvalues Ei and eigenvectors \Ei) of the Hamiltonian matrix were determined 
numerically using the appropriate subroutine for diagonalizing real symmetric 
matrices in the GNU Scientific Library [23] . 



4.2 Energy Eigenstates 

Fig. Q] shows a plot of the double well system V(x) — Vq(~2x 2 + x 4 ) with well 
depth Vo = 4.91345043. It also shows energy levels of the interacting system 
with the contact interaction chosen to be attractive (so the amplitude Uq is 
negative) . This can be achieved by tuning a homogeneous magnetic field to the 
Feshbach Resonance of the alkali metal atoms [24] [19] . 

Figs. [Ha through [He are plots of the ground state, \Ei), first excited state, 
\E 2 ), and third excited state \E±) of the interacting system, respectively. Each 
figure shows a plot of the probability density \(xi, X2\Ej)\ 2 , as well as cross 
sections of the wavefunctions (xi\Ej) for x 2 = [—1,0, 1]. The bosonic character 
of the states is evident from the fact that they are symmetric under x\ «-> x 2 
exchange. 

We can also compute the phase space distribution of the energy eigenstates 
and compare this with the classical surfaces of section in Fig. [2] . A phase space 
distribution of quantum states was first constructed by Wigner [25] . A smoothed 
version of the Wigner distribution was introduced by Husimi [26] and has proved 
particularly useful for comparison of classical and quantum phase space distri- 
butions. In the x — p representation, the Husimi Function Fh(xi,X2,Pi,p2)foT 
a symmetrized 2-particle wavefunction ipEj{x\,x 2 ) is defined as 

/ — — — — \ 1 [ dx x dx 2 
F h (x 1 ,x 2 ,pi,p 2 ) = / y Ej {xi,x 2 ) 

0\0 2 -K J V27T V27T 

(=°l-fl) 2 (x 2 -x 2 ) 2 
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where (x~i,p~i) and (x~2,p~2) are the centroids of the Gaussian wave packets in 
the phase space. 

In order to calculate the standard deviations <j\ and a 2 , we follow the 
same basic prescription as is normally followed for one-dimensional Husimi func- 
tions [28j , where the Gaussians appearing in the Husimi function are interpreted 
as harmonic oscillator coherent states. Therefore the standard deviation is the 
same as that of the harmonic oscillator, with the modification that the "stiff- 
ness" ujd be ^ where T is the double well period of motion, making cr, = ^J-^j^r 

(i is 1 or 2). This is generalized to two particles by choosing the single-particle 
energies that, when added up, come closest to the 2-particle energy. We then 
proceed to calculate u>d for each particle as a single noninteracting particle with 
the chosen energy level. A straightforward integration of the classical double 
well problem shows [UJ 

f^kife If ^<° (12) 

where K(k) is the complete elliptic integral of the first kind, k is given by 



^ = ^'\ + ^ , (13) 
and k 2 = 1 — ft 2 . 

Now that the standard deviations can be calculated, we have the prescrip- 
tion for numerically evaluating the full Husimi function. We cannot sketch the 
full four- dimensional function realistically, of course. However, we can sketch 
a "quantum Poincare map" of the Husimi function by strobing a particular 
value of X2 given p~2 > with the energy classically conserved at the quantum 
eigenvalue. Thus, we can plot 

fh{xi,pi) = F h (xi,x 2 = +l,_pi,p 2 = p2(xi,pi,Ej)), (14) 

where p 2 is determined from the condition at the unperturbed Hamiltonian 
H = Ej for a particular eigenstate of energy Ej . 

Figs. E]a through Oc are Husimi plots of the double well system for states 
I -El), \E 2 ) and \Ea). They can be compared with the corresponding classical 
Poincare sections in Figs. [Ha through Oc. The Husimi plot of \E\) is provided 
on Figs. [5]a.l and [5^.2. We notice that the highest probabilities are located 
in the separatrix region, where there is significant chaos in the classical map 
(Figs. HJa). However, there is a significant probability for the system to tunnel 
from the separatrix region to the interior near the well minima. This is where 
the bifurcated trajectories occur in Figs. dJa.l and[Ua.2. All the Husimis are 
symmetrical, apart from the interaction resonance in each case, under phase 
space inversion (x\— > — X\, p±— > — p\). The tunnelling probability from the 
separatrix into the wells is considerably reduced in the Husimi plot of \E±) 
(shown in Figs.Oc.l and[5]c.2). 
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5 Quantum Dynamics of the Driven System 



In order to control transitions between energy states of the two boson system, 
we drive the system with two sequential pulses of maser radiation with carrier 
frequency, u>f (u> s ) for the first (second) pulse. These frequencies are determined 
by the transitions of interest. The masers are projected in the dimension of 
particle confinement, thus the spacial dependence of the electric field 

Ei(x,t) = E 0i e^ kiX -^ + h.c (15) 

(i = f,s and j = y— 1) in that direction can be treated as linear, given that 
the wavelength of the radiation pulse(s) are in the microwave range, and the 
trapping length scales are ~ lOOnm in the double well. Ignoring purely temporal 
terms that only contribute an overall phase, the interaction Hamiltonian —Dm 
E(x, t) simplifies to jEoiDkiXe 3Wit + h.c, where D is the atomic dipole moment. 
Thus, the Hamiltonian of the driven system can be written as 

H(t) = H + [e f (t)sm(u} f t) + e 5 (i) sin(w s i)](>i + x 2 ), (16) 

where H is the Hamiltonian of the non-driven system in Eqn. ([7]) and the am- 
plitudes 6j(t) (i — f, s) of the driving fields have Gaussian shape. 

ei(t) = Ae^-^ 2 for i = f, s, (17) 

where Af (A s ) is the maximum amplitude of the first (second) pulse, and the 
dipole moment D and kts have been absorbed into the AiS. The duration of 
each pulse is controlled by the parameter (3 = , where r is a measure of the 
width of each pulse (similar to standard deviation of the Gaussian) . The time 
at which the maximum of the «th pulse occurs is U . 

The Schrodinger equation for the 2-boson system in the presence of the 
driving field is 

i^\m)=H{t)\m)- as) 

Given our numerical expressions for the energy eigenstates \Ej) of the non-driven 
system, we can expand \ip(t)) in terms of these states so that 

i 

The Schrodinger equation can then be written in the form 

= —iEjCj(t) — i[ef(t)sm(cjft) + e s (t) sin(ui s t)] ^ djj>Cj/(t), (20) 

f 

where Cj(t) — (Ej\ip(t)} is the probability amplitude to find the system in state 
\Ej) at time t and djji = {Ej\{x\+X2)\Eji) denotes dipole matrix elements taken 
with respect to the exact energy eigenstates of the undriven system. Values of 
the dipole matrix elements are given in Table [TJ 
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0.001 ... 
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Table 1: Dipolc Matrix elements for V a = 4.91345043 and U = -1.0. The first 
4X10 values are shown here. 



5.1 Floquet States 

For the case when the amplitude of the Gaussian pulses changes very slowly 
relative to the period of the carrier frequencies of the pulses, it is possible to 
use Floquet theory to study the dynamics of the driven system. As was shown 
in [13], [14], [16], one can divide the time over which the pulses act into a 
sequence of time intervals. During each time interval, the amplitude of the 
pulses is essentially constant while the carrier waves undergo many oscillations. 
Consider the time window centered at time t = tfi X . The Hamiltonian describing 
the dynamics during this time can be written 

H(t;t fix ) = H+ [e f (tf ix )shi{ujft) + e s (t fix ) sin(u; s t)] (xi +x 2 ). (21) 

If the two frequencies ujf and ui s are commensurate so that -f- = ^- where 
rif and n s are integers, then the Hamiltonian H(t;tfi X ) is time-periodic with a 
period 

T = tt(^- + — ). (22) 

LOf LO s 

Because H(t;tfi x ) is time periodic, Floquet theory can be used to analyze the 
dynamics of the system during the time window centered at t — tfi X . 

Let us assume that the Schrodinger equation, i-§i\ip(t)) = H(t;tfi X )\ip(t)) 
has a solution of the form 

m)) = e-^l^(t)), (23) 

where the state \4> a (t)} is time-periodic with period T and the phase fl a is 
real. If this is substituted into the Schrodinger equation we obtain the following 
eigenvalue equation 

{H(t;t fix ) - i^-)\4> a {t)) = fia|&»(*)>. (24) 



10 



The state \4> a {t)) is the ath Floquet eigenstate and Q a is the ath Floquet eigen- 
phase. The quantity HF(t)=H(t;tfi X )—i-j^ is a Hermitian operator and is called 
the Floquet Hamiltonian. The Floquet eigenstates \4> a (t)) form a complete or- 
thonormal basis and the Floquet eigenphases f2 a are conserved quantities [17] . 

The state of the boson system at time t can be expanded in a Floquet spectral 
decomposition as follows 

\m) = E^e-^'l^)) = ^(0 Q (O)|^(O))e-^ 4 |0 Q W). (25) 

a a 

Because the Floquet eigenstates are time-periodic, the state of the system at 
time t = T is given by 

\i>{T)) =^e- iQ « T |^(0))(^(0)^(0)). (26) 

a 

The Floquet evolution operator is therefore given by Uf(T) where 

U F (T)=J2^ maT \M0))(M0)\ (27) 

a 

and is a unitary operator. When the operator Uf(T) acts on the state of the 
boson system it evolves it forward in time by one period of the driving field. 

It is possible to compute the matrix elements of the Floquet evolution op- 
erator using energy eigenstates of the undriven system as the basis functions. 
Thus, 

Ujj.iT) = {E 3 \U F {T)\E r ) =Y,e- in * T (E j \<t> a (0))(cf> a (0)\E j '). (28) 

a 

The ath eigenvalue of this matrix is e _!n ° T and the ath eigenvector is given 
by a column matrix with entries (Ej\(j> a (0)). Since we will only have numerical 
expressions for the eigenvalues e _l0 " T , W e can only determine the value of the 
eigenphases fl a modulo 2ir. 

The numerical computation of the Floquet matrix Ujji (T) is achieved as fol- 
lows. Each column of the matrix can be constructed by solving the Schrodinger 
equation (with Hamiltonian H(t; t fix)) for one period T . Each column of the ini- 
tial state starts with a single entry (Ej\<p a (0)} = 1 for j = a and (Ej\tfi a (0)) = 
otherwise. The integration is done N times with a ranging from a = 1 to 
a = N. The numerical integrations were performed using the appropriate sub- 
routine for the 4" 1 order Runge-Kutta-Fehlberg method [22 from the GNU 
Scientific Library [23] . Each different initial condition yields one column of the 
Floquet matrix at time t = T. Performing these N integrations yields an TVx N 
Floquet evolution matrix. Numerically diagonalizing this matrix gives us the 
Floquet eigenphases and eigenstates. The numerical diagonalization of the non- 
Hcrmitian Floquet matrices were performed using the appropriate routine in 
the IBM™ Engineering and Scientific Subroutines Library (ESSL) [29]. This 
process can be performed for each value of t fi X and the resulting eigenphases 
and eigenstates plotted as functions of t fi x . 
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In order to determine the appropriate numerical truncation for the evaluation 
of the Floquet matrices, we used iteratively increasing values of N and checked 
the components of the Floquet states at the value(s) of tu x when the STIRAP 
amplitudes were the largest ie at t — tf or t s (see Eq. l[T7|l) until the higher 
components were too small to contribute to the dynamics. We chose TV = 25 as 
the final truncation. 

In the subsequent section we show that coherent transitions between sym- 
metrized two-particle boson states can be achieved for this system. Because of 
the sparsity on nonzero dipole matrix elements, the simplest transition, induced 
by the laser pulses, is from the ground state \Ei) to the fourth level l-Ej), via the 
intermediate state l-E^)- We show the behavior of the system for three different 
amplitudes of the radiation pulses. 

6 Case 1: STIRAP Ladder, First Pulse 2^4, 
Second Pulse 1 2 

Fig. [T] shows the energy levels of the double well system for well depth Vq — 
4.91345043 and interaction strength Uq — 1.0. The value of Vb was chosen so the 
radiation pulses would have carrier wave frequencies u>f and ui s commensurate 
with each other and so that ojf (oj s ) would be equal to the energy spacing E^—E2 
(Eq, — Ei), with a high degree of precision. The energy levels shown in the figure 
are the exact energy eigenvalues of the undriven two-boson symmetrized system. 

We plan to use radiation pulses to induce a coherent transition of the two 
boson system from its ground state \E\) to the excited state \E±). At t = the 
first pulse connects the levels E2 and £4 with zero detuning. The second pulse 
connects the levels E\ and E^. The ratio = I to eight decimal places. 

The dipole moments of these transitions have very different values. The 
dipole moment diA = {E 2 \{x\ + ^2) that couples the states I-E2) and \E&) is 
two orders of magnitude smaller than dipole moment = (Ex\(xi + X2)\E2) 
that couples the states |£i) and \E 2 ) (See Table 1). From Eqs. |16|) and (fTTj) . 
the amplitude of the first pulse is given by Afd 2 A and the amplitude of the 
second pulse is given by A s d\ t 2- Because the dipole coupling of the first pulse 
is so much smaller than that of the first pulse, we will make the electric field 
amplitude, Af, of the first pulse considerably larger than that of the second 
pulse, A s so that 

60^4^2,4 = 4^1,2- (29) 

The amplitudes for the two radiation pulses are plotted in Fig. [51 

The duration of each pulse can be controlled by varying the pulse width 
parameter r. We let ttot denote the total time over which both pulses act on 
the system. We choose the following values for the pulse parameters 

11 2 

t = -hot, tf = -Uot, and t s = -t to t- (30) 
00 o 

In the sections below, we will study the effect of these radiation pulses on 
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the boson system for two values of eo- In both of these cases, we set a value for 
eo and set a suitable truncation value N for the Floquet evolution matrix. 

The Floquet eigenphases lie within a fundamental zone (they are determined 
modulo lu) that is taken to be {0,w} where ui = ^- — 1.05303. They can 
be plotted as a function of tf ix . For closely spaced values of tf ix , Floquet 
eigenstates belonging to different eigenphases, at neighboring values of will 
be orthogonal. This can be exploited to tag and follow the evolution of each 
eigenstate and eigenphase as a function of t fix ■ 

In subsequent sections, we label each Floquet eigenphase based on it's dom- 
inant dependence on the undriven Hamiltonian eigenstates, \Ej) at tfi x = 0. 
For the three levels, \Ei), \E 2 ) and \E±), that are connected by the STIRAP 
pulses, the corresponding Floquet eigenstates have the following structure and 
labels: 

1. The eigenphase whose corresponding Floquet eigenstate is dominated by 
the undriven ground state \Ei) at tfi X — is labelled as and the 
Floquet eigenstate as \4>a)- 

2. The eigenphase whose corresponding Floquet eigenstate is dominated by 
the undriven state [I-E4) — I-E2)] at tf ix — is labelled as Us and the 
Floquet eigenstate as |^b). 

3. The eigenphase whose corresponding Floquet eigenstate is dominated by 
the undriven state ^= [\E±) + I-E2)] at tf ix — is labelled as fic and the 
Floquet eigenstate as \<j>c)- 

4. The eigenphase whose corresponding Floquet eigenstate is dominated by 
the undriven state \Ej) at tji X = is labelled as SId and the Floquet 
eigenstate as \4>d)- 

The symmetric and antisymmetric Floquet states [\E^)±\E2)] are induced 
by the first radiation pulse (which couples the states \E 2 ) and I-E4)) even though 
the amplitude of the first radiation pulse may be very small. 

The results of the Floquet analysis described above can be compared to the 
exact dynamics of the system obtained by solving the full Schrodingcr equation 
in Eq. (j2"0)) for the exact state of the driven system \t/j(t)). In solving for \ip(t)}, 
we will always start at time t — with the system in the ground state |^(0)) = 
\E\) of the undriven Hamiltonian. We can then plot the probability Pj{t) = 
\{Ej\^{t))\ 2 of finding the syste m in the undriven energy level \Ej) as a function 
of time t for various values of ttot- We cannot show strobed Husimi plots of the 
Floquet states as they evolve across tfi X , nor can we show classical Poincare 
maps of the system during those times, since the system has five degrees of 
freedom during the STIRAP process. 

6.1 Pulse Amplitude eo = 10 

Fig. 0a we plot the Floquet eigenphases as a function of tfix in units of the total 
pulse time ttot- The eigenphases of interest are the ones involved in the STIRAP 
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process (ie QaMbMc) which lie in the interval {0.488^0.5} in Fig. 0a. Fig. 0b 
shows a magnification of that region. It is clear that the three levels contribute in 
a manner characteristic of a traditional STIRAP ladder process approximated by 
a three-level system [TJ]. A three-level avoided crossing occurs at tfi x ~ 0.5i tot , 
and a coherent population transfer takes place from the ground state to the third 
excited state. This is further confirmed by Figs. [5]a and [5]b. Fig. [5]a shows 
the evolution of the dependence of \4>a) on the undriven energy eigenstates \Ej) 
as a function of time t/j x - A population transfer occurs from the ground state 
(labelled " 1" ) to the fourth energy level E 4 (labelled "4" ) . 

Fig. [5]b shows the actual time evolution of the state of the system \ip(t)) 
obtained by solving the Schrodinger equation as a function oft for i tot = 28000. 
In Fig.[8lb, we plot the value of Pj{t) = \(Ej\ip(t))\ 2 as a function of time t. The 
real time evolution is very close to the evolution of the Floquet eigenstate \4>a) 
as a function of tfi x . This indicates that the evolution is governed by a single 
Floquet eigenstate and that the process is adiabatic. The small oscillations and 
deviations of Fig. Ob from Fig.JSja can be attributed to nonadiabatic effects [50] . 

6.2 Pulse Amplitudes e = 115 

We now set eo at a higher value of 115. Fig. [!|]a shows the evolution of the 
Floquet eigenphases as a function of tfi x . Fig. Ob shows a magnification of the 
region containing eigenphases fi^, f2#, £lc and Qjj. We notice the prominence 
of a new Floquet state \4>d) (with corresponding eigenphase profile Flo) which, 
at tfi x — 0, is displaced in value from f^, fig and Qq- At tfi X = 0, \4>d) is 
dominated by \Ey), This occurs due to the near-resonance between the 2 — 4 
transition and the 4 — 7 transition (see Fig. [l}. 

We also note, from TableQ] that the 4 — 7 dipole moment (0.015) is an order 
of magnitude higher than the dipole moment of the connected states 2 — 4. 
Since the 2 — 4 resonance is very close to the 4 — 7 resonance, the evolution 
of the corresponding eigenphase fiu is affected and influences the evolution of 
the eigenphase Ha- However no measureable avoided crossing with £Ia occurs. 
Therefore, Qa and £Id appear to cross (or under go an avoided crossing very 
closely spaced) and do not contribute anything significant to the dynamics. 

There are three avoided crossings that can affect transitions in the system 
for eo = H5. First, Qa appears to undergo an avoided crossing with £lc at 
tfi x — 0.29ti O t (see Fig. 10. a. 1). Then the same pair of states avoid each other 
again at tfi X ~ 0.35tt o t (see Fig. 10. b). Finally, the states FIa, ^b, and flc 
undergo the standard STIRAP transition at tfi x ~ 0.5tt o t- The dependence of 
\</>a),\4>b), and \4>c) on the unperturbed energy eigenstates is shown in Fig. [TT1 
The influence of these avoided crossings on these states is clearly seen. These 
avoided crossings are manifestations of classical chaos in the quantum dynamics 
as elaborated in the introduction. 

We can now use known properties of avoided crossings to analyze this process 
in more detail. When two Floquet eigenphases Q a and £lp approach and undergo 
an isolated avoided crossing, the probability P a p that the system switches from 
one Floquet state to the other can be calculated from the Landau-Zener formula 
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for two-level systems [3T] (note that use of this estimate for multi-level systems 
assumes that other levels are not significantly involved in the avoided crossing) . 
In our dimensionless units, the Landau-Zener probability P a $ is 



P a p = exp 



2T 



a/3 



(31) 



where 5Q a p is the (minimum) spacing between £l a and fig at the avoided cross- 
ing and T a p is the magnitude of the rate of change (slope) of the Floquet 
eigenphases in the immediate neighborhood of the avoided crossing. Thus, 



a/3 



<Kl a d£l 



dt dt 



(32) 



where is the slope of the eigenphase curve fl a in the neighborhood of the 
avoided crossing. If the system switches between the two Floquet states at the 
avoided crossing, (if P Q( g«l), then the energy eigenstates of the undriven system 
which contribute the evolution do not change significantly. On the other hand, 
if P a p~0, then the system follows a single Floquet state through the avoided 
crossing, but there can be significant change in the energy eigenstates of the 
undriven system that contribute to the dynamics. 

The value of T a p depends on the duration of the pulses ttot because that 
determines the slopes of the Floquet eigenphase curves as they enter and leave 
the avoided crossing. To make this explicit, we can write 

dn 1 dn 

— = — , (33) 

dt ttot dT ac 

where r ac = tfi X /t to t is the time (normalized to t to t) at which the avoided 
crossing occurs. The quantity f = 4j^- has very weak dependence on t tot . We 
can now write the Landau-Zener transfer probability in the form 

P a/3 = exp [-ttot1a,p] ■ (34) 

The quantity 7 Qj/ 3 = T ^|!°^ has weak dependence on t tot . Thus, the transfer 

probability P a p will be very small if t to t > l/7a,/3- 

We can compute the Landau-Zener probability for the avoided crossings that 
occur for this case. We use the information from Fig. [10] to calculate 7 a „8- For 
the avoided crossing between f^c and £Ia, shown in Fig. [Tola at tfi X ~ 0.29, 
Sn c ,A = 0.000605 and Y c ,a = 0.0929395. Therefore, j c ,a = 6.19xl0~ 6 and 
we must have ttot > 1.62xl0 5 to have a small probability that the system will 
transfer from one Floquet state to the other. For the second avoided crossing 
between ile an d &a shown in Fig. [TOlb at tf i X ~ 0.35, 5flc\A — 0.0032651 and 
Tc,A = 0.136394. This means that t to t > 8144.8 for an adiabatic passage (no 
change in Floquet eigenstate) through the avoided crossing. 

Fig.[12]a, which has a relatively small value of ttot (i*oi~700), shows no effect 
of these first two avoided crossings, although it does show the effect of the three- 
way avoided crossing that occurs about halfway into the total time. For this 
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case, the pulses appear to leave the system in a superposition of Floquet states 
\4>a) and \4>c)- The effect of the avoided crossing at i/ i:E ~0.29 is also absent 
in Fig. rj2] b, where t to t = 7000. However, the effect of the avoided crossing at 
tfix ~ 0.35 can be seen in the Fig. A complex mixing of states \Ej) occurs 
just after this avoided crossing as the large central three-way avoided crossing 
comes into play, and in the end the system is again left in a superposition of 
Floquet states \4>a) and \(j>c)- The avoided crossing at £/; x ^0.29 finally starts 
to manifest itself in Fig. [T2] (c), where ttot — 72,000 and is clearly visible in 
Fig.[TJ]d, where t to t = 720, 000. Indeed, in Fig.[T2]d, the system follows a single 
Floquet state through the entire process. This is confirmed by comparing the 
evolution of the Floquet eigenstate \4>a) in Fig. [TT]b to the exact time evolution 
in Fig. [TJJd. They are essentially identical. 

7 Case 2: STIRAP Ladder, First Pulse 2^4, 
Second Pulse 1 <-> 2 and nearly tuned to 4 <-> 7 

We now want to show an interesting effect that can occur in a multilevel system. 
We adjust the shape of the double-well potential so that there is a resonance 
between the 2^4 and 4^7 transitions that is almost exact (to within 10~ 3 
units of energy). Fig. [13] shows the energy levels of the double well system 
for wells that are a little deeper than in Case 1. Here, Vq = 7.2912229 and 
interaction strength Uq = —1.0. The energy levels shown in the figure are the 
exact energy eigenvalues of the undriven two-boson symmetrized system. 

The classical dynamics of the system, for these deeper potential wells, is 
qualitatively the same as in Case 1, depicted in Figs. [2] and Figs. [31 However, 
deepening the wells lowers the quantum energies. Nonetheless, the classical 
dynamics that was observed at energies E\ , E2 and £4 in Case 1 is very similar 
to that seen for the corresponding energies Ei, E 2 and E4 for Case 2. For Case 
2, that we consider in this section, the energy Ej will play a significant role. 
In Fig. [Mja.l { x i < 0) and [T4la.2 (xi > 0) we show surfaces of section for 
the classical non-driven interacting system for an energy equal to the seventh 
energy level E? of Case 2. The chaos is more spread out and the effect of a new 
bifurcation can be seen. 

We use radiation pulses to induce a coherent transition of the two-boson 
system from its ground state in a manner similar to Case 1. At t = the 
first pulse connects the levels E 2 and E4 with zero detuning. The second pulse 
connects the levels E\ and En. In this case, the ratio — = As in Case 1, 
there is a significant difference (3 orders of magnitude) between dipole moments 
di2 and d 2 4 (see tabled]). Therefore, the peak amplitudes of the first and second 
pulses are adjusted in accordance with Eq. l[2"9"]) (see Fig. [5]). In Fig.rjSJa we plot 
the Floquet eigenphases as a function oltf ix in units of the total pulse time i tot . 
The fundamental zone has been set to {— ~, j}, where w, the commensurate 
frequency, is given by -=r and T is calculated from Eq. ([22]) . In this case, the 
ratio rf- = The eigenphases of interest are labelled in the same manner 
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0.0002 




0.0004 





0.002 




Table 2: Dipole Matrix elements for V = 7.2912229 and U = -1.0. The first 
4X10 values are shown here. 



as in case 1 ie as VIa, ^b,^c, and fin- These eigenphases lie in the interval 
{—0.02, +0.01} in Fig. [T^la. Fig. Qj}]b shows a magnification of that region. 

As we can see from Fig. rT5]b, all four eigenphases CIa, &b, and Qd, 
participate in a complicated set of avoided crossings. Firstly, an avoided crossing 
between flc and £Id very near tu x = causes them to switch their supports. 
Thus, \4>c) is predominantly supported by \E~j) after this crossing is avoided. 
The next avoided crossing of importance is the one between ft a and flc at 
tj ix ~ 0.50t tot , which causes the support of \(/>a) to change from \E{) to that 
of \4>c) viz- \Et). Thus, a complete population transfer from \E{) to \Ei) is 
possible. The effect of these avoided crossings and the possible behavior of 
the system, as the radiations pulses pass through the system, can be seen in 
Fig. HH Figs. HUa through [TBld show the time strobed plots of \(/)a} through 
|</>d), respectively, analogous to Fig. QT] of Case 1. It is clear from Fig. [TBl a 
that this unexpected transition from energy level E\ to E-j should be possible 
to achieve, producing a marked influence of classical chaos in the quantum 
dynamics as elaborated in the introduction. 

In order to obtain a rough estimate of the pulse time t tot needed to achieve 
true adiabatic behavior, we apply the Landau Zener formula Eqn. (|3~4"|) in the 
same manner as previously done in Case 1, even though it was not meant to 
be applicable when multiple avoided crossings are involved. For the avoided 
crossing of VIa with Qc at tf ix ~ 0.50t to t, shown in Fig.[l5]b, SVIa^c — 0.0047117 
and f A ,C = 0.0403021. Therefore ~/ A ,c = 0.000865261 and we must have t tot > 
1.2xl0 3 to have a small probability that the system transfers from \4>a) to \4>c)- 

Fig.[T71a through fTTlc show the actual time evolution of the system, starting 
from the ground state \E\) at t = for different values of ttot- For small values of 
t to t, below the threshold calculated with the Landau Zener formula (Fig. [T7la). 
we see a partial coherent transfer to I-E4), as demonstrated above. When ttot 
is well above threshold, as in Fig.[l7lc, complete population transfer to \E-f) is 
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achieved and we appear to have reached approximately adiabatic behavior. 

8 Conclusions 

We have observed some unusual behavior in the classical and quantum dynamics 
of two bosons in a double well. Chaos in the separatrix region of the classical 
version of the coupled system corresponds with regions of high probability in 
the quantum Poincare map. However, a noticeable tunnelling has been observed 
from the separatrix into the individual wells. We have also demonstrated the 
feasibility of a controlled excitation of the system into a higher energy state using 
STIRAP. The STIRAP pulses destroy symmetries and produce chaos that we 
can detect by observing avoided crossings in the Floquet eigenphase spectrum. 
The chaos produced by additional resonances produce avoided crossings that 
can cause coherent population transfer to higher states. Thus, radiation pulses 
can be used to exert coherent control of the coupled boson system through 
chaos assisted adiabatic passages, just as has been recorded for systems with 
lower degrees of freedom. 
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Figure 1: Plot of the double- well potential experienced by each boson in case 1. 
All units are dimensionless. The energy levels, E\ = —3.71958, E 2 — —2.66655 
and £4 = 2.5986 of the interacting two-boson system (interaction strength Uq = 
— 1.0) are also sketched, with wavy arrows denoting the levels connected by the 
STIRAP pulses. Here, V a = 4.91345043. 



Figure 2: Classical Poincare maps (of (pi,x\) for x 2 = 1.0 and p 2 > 0) for two 
interacting particles in the double well potential. All units are dimensionless. 
Here Vq = 4.91345043 and Uq = —1.0. The interaction is approximated by 
an attractive Gaussian potential of width a c = 0.05. (a) Energy E = E\ — 
-3.71958. (b) Energy E = E 2 = -2.66655. and (c) Energy E = E t = 2.5986. 
A unit area of the phase space equals h 



Figure 3: Classical Poincare maps for lower energies with the bifurcating res- 
onance magnified. All units are dimensionless. A unit area of the phase space 
equals ft. (a) Energy E = E 1 = -3.71958. (b) Energy E = -3.8, (c) En- 
ergy E = -4.0, (d) Energy E = -4.5, (e) Energy E = -4.7, and (f) Energy 
E = —5.0. Note the increased prominence of the smaller resonance as the en- 
ergy decreases from E\ , as well as the bifurcation in the other resonance as the 
energy increases from (e). 



Figure 4: Plots of energy eigenfunctions for the two interacting bosons in a 
double well potential. All units are dimensionless. (al) Contour plot of the 
probability density |(xi,a;2|£a)| 2 • (a2) The cross-section of the wavefunction 
at x 2 = —1.0. (a3) The cross-section of the wavefunction at x 2 = 0. (a4) 
The cross-section of the wavefunction at x 2 = +1. (bl) Contour plot of the 
probability density \(xi,x 2 \E 2 )\ 2 . (b2) The cross-section of the wavefunction 
at x 2 = —1.0. (b3) The cross-section of the wavefunction at x 2 = 0. (b4) 
The cross-section of the wavefunction at x 2 = +1. (cl) Contour plot of the 
probability density \(x\, x 2 \E4)\ 2 . (c2) The cross-section of the wavefunction at 
x 2 = —1.0. (c3) The cross-section of the wavefunction at x 2 = 0. (c4) The 
cross-section of the wavefunction at x 2 = +1. 



Figure 5: Husimi functions for several energy eigenfunctions for two interacting 
bosons in a double well potential, with Vb = 4.91345043 and Uq = —1.0. All 
units are dimensionless. The figures show density plots of the Husimi distribu- 
tion in the (xi,pi) plane fixed at x 2 = 1.0 and p 2 > and subject to classical 
energy conservation, (a) The energy eigenstate \E\). (b) The energy eigenstate 
\E 2 ). (c) The energy eigenstate \E 4 ). 
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Figure 6: The STIRAP pulse amplitudes as a function of time for case 1. The 
first pulse (in time) connects the intermediate state to the final state of the 
STIRAP process. All units are dimcnsionlcss. The second pulse (in time) 
connects the initial state and the intermediate state. The total time t to t is chosen 
arbitrarily, but the centroids of the pulses are kept at = = \ijf~ t = \ 



Figure 7: Floquet eigenphase plots for case 1 with eo = 10 and u> = 1.05303. 
All units are dimensionless. (a) A plot of the Floquet eigenphases as a function 
of t f ix /t to t- (b) Magnification of the region of interest of in (a). 



Figure 8: (a) Plot of |(E^|^>a)| 2 as a function of tfi X /t to t for case 1 with e = 
10. All units are dimensionless. The numbers attached to each curve indicate 
the particular eigenstate \Ej) represented, (b) The exact time evolution of 
\(Ej\%p(t)}\\ 2 , obtained by solving the Schrodingerer equation with initial state 
IV'(O)) = The total pulse time is t tot = 24000. 



Figure 9: Floquet eigenphase plots for case 1 with e = 115 and w = 1.05303. 
All units are dimensionless. (a) A plot of the Floquet eigenphases as a function 
of t fi X j ttot- (b) Magnification of the region of interest of in (a). 



Figure 10: Magnification of avoided crossings in Fig. 12(a). All units are di- 
mensionless. (a) The first avoided crossing is between the curves for eigenphases 
Qa and Clc- (b) A second avoided crossing between and f^c that appears 
to reverse the effects of the first avoided crossing. 



Figure 11: (a) Plot of |(£^|<^c*)| 2 as a function of tfi X /t to t for case 1 with 
eo = 115. All units are dimensionless. The numbers attached to each curve 
indicate the particular eigenstate \Ej) represented, (b) Plot of K-EjI^a)! 2 as a 
function oitfi X /t to t- (c) Plot of \(Ej\(f>B)\ 2 as a function of tfi X /t to t- 



Figure 12: Plots of \(Ej\tp(t)}\\ 2 as a function of time for case 1 with e = 115 
and initial condition |i/j(0)}| = for different values of t to t- All units are 

dimensionless. (a) t tot = 72. (b) t tot = 7000. (c) t tot = 24000. (d) t tot = 72000. 
(e) ttot = 720000. 
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Figure 13: Plot of the double- well potential experienced by each boson in case 2. 
All units are dimensionless. The energy levels, E\ = —6.42262, E% = —5.68883 
and Ei — 0.640055 of the interacting two-boson system (interaction strength 
Uq = —1.0) are also sketched, with wavy arrows denoting the levels connected 
by the STIRAP pulses. Note the slightly detuned resonance between the 2^4 
and the 4 <-> 7 levels where E 7 = 6.96998. Here, V = 7.2912229. 



Figure 14: Classical Poincare maps (of (pi,x\) for x 2 = 1.0 and p 2 > 0) for two 
interacting particles in the double well potential. All units are dimensionless. 
Here Vq = 7.2912229 and Uo = —1.0. The interaction is approximated by an 
attractive Gaussian potential of width a c = 0.005. Energy E = E-j = 6.96998. 
A unit area of the phase space equals Ti. 



Figure 15: Floquet eigenphasc plots for case 2 with 6q = 115 and lo = 
0.091722994. All units are dimensionless. The fundamental zone here is 
{~^T' f } ( a ) -^ plot of the Floquet eigenphases as a function of tf ix /t to t- (b) 
Magnification of the region of interest of in (a) . 



Figure 16: (a) Plot of K-Ejl^)! 2 as a function of tfi x /t to t for case 2 with 
e = 115. All units are dimensionless. The numbers attached to each curve 
indicate the particular eigenstate \Ej) represented, (b) Plot of \{Ej\4>B)\ 2 as a 
function of tfi X /t to t- (c) Plot of \{Ej\4>c)\ 2 as a function of tfi X /t to t-{&) Plot of 
\(Ej\(j) D )\ 2 as a function of t fix /t tot . 



Figure 17: Plots of |(.Ey|^;(£))|| 2 as a function of time for case 2 with e — 115 
and initial condition ^(O))! = \Ej)\, for different values of t tot . All units are 
dimensionless. (a) t tot = 1000(b) t tot = 3600. (c) t tot = 36000. 
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